import pandas as pd
import sklearn
import numpy as np
import matplotlib.pyplot as plt
import os
import subprocess
import nibabel as nib
import nilearn.plotting as plotting
from nilearn.image import mean_img
import glob
from tensorflow import keras
import tensorflow as tf
from sklearn.model_selection import train_test_split
from sklearn.model_selection import cross_val_score
from sklearn.model_selection import cross_validate
from sklearn.model_selection import KFold
# Set FSL environment variables
os.environ['FSLDIR'] = '/Users/asaraog/fsl'
os.environ['FSLOUTPUTTYPE'] = 'NIFTI_GZ'
os.environ['PATH'] += os.pathsep + os.path.join(os.environ['FSLDIR'], 'bin')
# Define file paths and parameters
data_dir = 'ds002185-download'
2024-08-27 16:37:28.489095: I tensorflow/core/platform/cpu_feature_guard.cc:182] This TensorFlow binary is optimized to use available CPU instructions in performance-critical operations. To enable the following instructions: AVX2 FMA, in other operations, rebuild TensorFlow with the appropriate compiler flags.
#Example for EDA
subject_id = '01'
run_id = '01'
func_file = os.path.join(data_dir, f'sub-{subject_id}', f'func/sub-{subject_id}_task-odors_run-{run_id}_bold.nii.gz')
output_dir = 'preproc'
os.makedirs(output_dir, exist_ok=True)
#View image before preprocessing
img = nib.load(func_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='Raw')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order)
# Step 1: Motion Correction with MCFLIRT
motion_corrected_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_mcf.nii.gz')
mcflirt_cmd = f"mcflirt -in {func_file} -out {motion_corrected_file}"
subprocess.run(mcflirt_cmd, shell=True, check=True)
img = nib.load(motion_corrected_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='Motion Corrected')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order)
# Step 2: Extract a single volume from the 4D data for BET
extracted_volume_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_extracted.nii.gz')
extract_volume_cmd = f"fslroi {motion_corrected_file} {extracted_volume_file} 0 1"
subprocess.run(extract_volume_cmd, shell=True, check=True)
# Step 3: Brain Extraction
bet_output_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_brain.nii.gz')
bet_mask_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_brain_mask.nii.gz')
bet_cmd = f"bet {extracted_volume_file} {bet_output_file} -m -g 0.65"
subprocess.run(bet_cmd, shell=True, check=True)
# Step 4: Apply brain mask to the full 4D fMRI data
masked_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_masked.nii.gz')
apply_mask_cmd = f"fslmaths {motion_corrected_file} -mas {bet_mask_file} {masked_file}"
subprocess.run(apply_mask_cmd, shell=True, check=True)
img = nib.load(masked_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='After BET')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order) /usr/local/lib/python3.11/site-packages/nilearn/plotting/html_document.py:102: UserWarning: It seems you have created more than 10 nilearn views. As each view uses dozens of megabytes of RAM, you might want to delete some of them. warnings.warn(
# Step 5: Spatial Smoothing
smoothed_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_smooth.nii.gz')
smooth_cmd = f"fslmaths {masked_file} -s 2.548 {smoothed_file}"
subprocess.run(smooth_cmd, shell=True, check=True)
img = nib.load(smoothed_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='Smoothed')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order) /usr/local/lib/python3.11/site-packages/nilearn/plotting/html_document.py:102: UserWarning: It seems you have created more than 10 nilearn views. As each view uses dozens of megabytes of RAM, you might want to delete some of them. warnings.warn(
# Step 6: High-Pass Temporal Filtering
highpassed_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_highpass.nii.gz')
highpass_cmd = f"fslmaths {smoothed_file} -bptf {125 / 2.3548} -1 {highpassed_file}"
subprocess.run(highpass_cmd, shell=True, check=True)
img = nib.load(highpassed_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='High-Temporal Pass')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order) /usr/local/lib/python3.11/site-packages/nilearn/plotting/html_document.py:102: UserWarning: It seems you have created more than 10 nilearn views. As each view uses dozens of megabytes of RAM, you might want to delete some of them. warnings.warn(
# Step 7: Compute Mean Image for Nonlinear Registration to MNI Template Space
mean_image_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_mean.nii.gz')
mean_cmd = f"fslmaths {highpassed_file} -Tmean {mean_image_file}"
subprocess.run(mean_cmd, shell=True, check=True)
# Step 8: Resample MNI Template to Match Data FOV
mni_template = '/Users/asaraog/fsl/data/standard/MNI152_T1_2mm_brain.nii.gz'
resampled_template_file = os.path.join(output_dir, f'mni_template_resampled.nii.gz')
resample_cmd = f"flirt -in {mni_template} -ref {mean_image_file} -out {resampled_template_file} -applyisoxfm 2.0"
subprocess.run(resample_cmd, shell=True, check=True)
# Step 9: Nonlinear Registration to MNI Template Space
mni_registered_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_mni.nii.gz')
nonlinear_registration_cmd = f"fnirt --in={mean_image_file} --ref={resampled_template_file} --iout={mni_registered_file} --fout={output_dir}/warp_field --jout={output_dir}/jacobian"
subprocess.run(nonlinear_registration_cmd, shell=True, check=True)
# Step 10: Apply the MNI Registration Transform to All Time Points
mni_registered_all_timepoints_file = os.path.join(output_dir, f'sub-{subject_id}_task-odors_run-{run_id}_bold_preprocessed.nii.gz')
apply_transform_cmd = f"applywarp --ref={resampled_template_file} --in={highpassed_file} --warp={output_dir}/warp_field --out={mni_registered_all_timepoints_file}"
subprocess.run(apply_transform_cmd, shell=True, check=True)
img = nib.load(mni_registered_all_timepoints_file)
data = img.get_fdata()
if data.ndim == 4:
img = nib.Nifti1Image(data[:, :, :, data.shape[3] // 2], img.affine)
plotting.view_img(img, threshold=None, title='MNI Registered All Timepoints')
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order) /usr/local/lib/python3.11/site-packages/nilearn/plotting/html_document.py:102: UserWarning: It seems you have created more than 10 nilearn views. As each view uses dozens of megabytes of RAM, you might want to delete some of them. warnings.warn(
import os
import subprocess
# Define base directories
data_dir = 'ds002185-download' # Input and output directory
# Get a list of all subject IDs
subject_ids = [d for d in os.listdir(data_dir) if d.startswith('sub-')]
for subject_id in subject_ids:
if subject_id == 'sub-21':
print(f"Skipping {subject_id} due to missing functional files.")
continue
subject_dir = os.path.join(data_dir, subject_id, 'func')
if not os.path.exists(subject_dir):
print(f"Skipping {subject_id} due to missing 'func' directory.")
continue
func_files = [f for f in os.listdir(subject_dir) if f.endswith('_bold.nii.gz')]
for func_file in func_files:
# Extract run_id from the file name (assuming format sub-01_task-odors_run-01_bold.nii.gz)
run_id = func_file.split('_')[2].split('-')[1] # Extract '01' or '02'
if subject_id == 'sub-22' and run_id == '01':
print(f"Skipping {subject_id} run-{run_id} due to corrupt file.")
continue
func_file_path = os.path.join(subject_dir, func_file)
output_dir = subject_dir # Output in the same place as input
os.makedirs(output_dir, exist_ok=True)
# Step 1: Motion Correction
motion_corrected_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_mcf.nii.gz')
mcflirt_cmd = f"mcflirt -in {func_file_path} -out {motion_corrected_file}"
subprocess.run(mcflirt_cmd, shell=True, check=True)
# Step 2: Extract a single volume from the 4D data for BET
extracted_volume_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_extracted.nii.gz')
extract_volume_cmd = f"fslroi {motion_corrected_file} {extracted_volume_file} 0 1"
subprocess.run(extract_volume_cmd, shell=True, check=True)
# Step 3: Brain Extraction
bet_output_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_brain.nii.gz')
bet_mask_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_brain_mask.nii.gz')
bet_cmd = f"bet {extracted_volume_file} {bet_output_file} -m -g 0.65"
subprocess.run(bet_cmd, shell=True, check=True)
# Step 4: Apply brain mask to the full 4D fMRI data
masked_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_masked.nii.gz')
apply_mask_cmd = f"fslmaths {motion_corrected_file} -mas {bet_mask_file} {masked_file}"
subprocess.run(apply_mask_cmd, shell=True, check=True)
# Step 5: Spatial Smoothing
smoothed_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_smooth.nii.gz')
smooth_cmd = f"fslmaths {masked_file} -s 2.548 {smoothed_file}"
subprocess.run(smooth_cmd, shell=True, check=True)
# Step 6: High-Pass Temporal Filtering
highpassed_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_highpass.nii.gz')
highpass_cmd = f"fslmaths {smoothed_file} -bptf {125 / 2.3548} -1 {highpassed_file}"
subprocess.run(highpass_cmd, shell=True, check=True)
# Step 7: Compute Mean Image for Nonlinear Registration to MNI Template Space
mean_image_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_mean.nii.gz')
mean_cmd = f"fslmaths {highpassed_file} -Tmean {mean_image_file}"
subprocess.run(mean_cmd, shell=True, check=True)
# Step 8: Resample MNI Template to Match Data FOV
mni_template = '/Users/asaraog/fsl/data/standard/MNI152_T1_2mm_brain.nii.gz'
resampled_template_file = os.path.join(output_dir, f'mni_template_resampled.nii.gz')
resample_cmd = f"flirt -in {mni_template} -ref {mean_image_file} -out {resampled_template_file} -applyisoxfm 2.0"
subprocess.run(resample_cmd, shell=True, check=True)
# Step 9: Nonlinear Registration to MNI Template Space
mni_registered_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_mni.nii.gz')
nonlinear_registration_cmd = f"fnirt --in={mean_image_file} --ref={resampled_template_file} --iout={mni_registered_file} --fout={output_dir}/warp_field --jout={output_dir}/jacobian"
subprocess.run(nonlinear_registration_cmd, shell=True, check=True)
# Step 10: Apply the MNI Registration Transform to All Time Points
mni_registered_all_timepoints_file = os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_preprocessed.nii.gz')
apply_transform_cmd = f"applywarp --ref={resampled_template_file} --in={highpassed_file} --warp={output_dir}/warp_field --out={mni_registered_all_timepoints_file}"
subprocess.run(apply_transform_cmd, shell=True, check=True)
# Clean up intermediate files
intermediate_files = [
motion_corrected_file,
extracted_volume_file,
bet_output_file,
bet_mask_file,
masked_file,
smoothed_file,
highpassed_file,
mean_image_file,
resampled_template_file,
mni_registered_file,
os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_mean_to_mni_template_resampled.log'),
os.path.join(output_dir, f'{subject_id}_task-odors_run-{run_id}_bold_mean_warpcoef.nii.gz'),
os.path.join(output_dir, 'warp_field.nii.gz'),
os.path.join(output_dir, 'jacobian.nii.gz')
]
for file in intermediate_files:
if os.path.isfile(file):
os.remove(file)
print(f"Preprocessing complete for {subject_id}_run-{run_id}! Intermediate files deleted.")
Preprocessing complete for sub-13_run-01! Intermediate files deleted. Preprocessing complete for sub-13_run-02! Intermediate files deleted. Preprocessing complete for sub-14_run-01! Intermediate files deleted. Preprocessing complete for sub-14_run-02! Intermediate files deleted. Skipping sub-22 run-01 due to corrupt file. Preprocessing complete for sub-22_run-02! Intermediate files deleted. Preprocessing complete for sub-15_run-02! Intermediate files deleted. Preprocessing complete for sub-15_run-01! Intermediate files deleted. Preprocessing complete for sub-12_run-02! Intermediate files deleted. Preprocessing complete for sub-12_run-01! Intermediate files deleted. Preprocessing complete for sub-08_run-01! Intermediate files deleted. Preprocessing complete for sub-08_run-02! Intermediate files deleted. Preprocessing complete for sub-01_run-02! Intermediate files deleted. Preprocessing complete for sub-01_run-01! Intermediate files deleted. Preprocessing complete for sub-06_run-02! Intermediate files deleted. Preprocessing complete for sub-06_run-01! Intermediate files deleted. Preprocessing complete for sub-07_run-01! Intermediate files deleted. Preprocessing complete for sub-07_run-02! Intermediate files deleted. Preprocessing complete for sub-09_run-02! Intermediate files deleted. Preprocessing complete for sub-09_run-01! Intermediate files deleted. Preprocessing complete for sub-17_run-02! Intermediate files deleted. Preprocessing complete for sub-17_run-01! Intermediate files deleted. Preprocessing complete for sub-10_run-02! Intermediate files deleted. Preprocessing complete for sub-10_run-01! Intermediate files deleted. Preprocessing complete for sub-19_run-01! Intermediate files deleted. Preprocessing complete for sub-19_run-02! Intermediate files deleted. Skipping sub-21 due to missing functional files. Preprocessing complete for sub-20_run-01! Intermediate files deleted. Preprocessing complete for sub-20_run-02! Intermediate files deleted. Preprocessing complete for sub-18_run-02! Intermediate files deleted. Preprocessing complete for sub-18_run-01! Intermediate files deleted. Preprocessing complete for sub-11_run-01! Intermediate files deleted. Preprocessing complete for sub-11_run-02! Intermediate files deleted. Preprocessing complete for sub-16_run-01! Intermediate files deleted. Preprocessing complete for sub-16_run-02! Intermediate files deleted. Preprocessing complete for sub-05_run-01! Intermediate files deleted. Preprocessing complete for sub-05_run-02! Intermediate files deleted. Preprocessing complete for sub-02_run-01! Intermediate files deleted. Preprocessing complete for sub-02_run-02! Intermediate files deleted. Preprocessing complete for sub-03_run-02! Intermediate files deleted. Preprocessing complete for sub-03_run-01! Intermediate files deleted. Preprocessing complete for sub-04_run-02! Intermediate files deleted. Preprocessing complete for sub-04_run-01! Intermediate files deleted.
#Example trial_type for one subject and the time to DICOM volume in the matrix
events_df = pd.read_csv('ds002185-download/sub-05/func/sub-05_task-odors_run-02_events.tsv', delimiter='\t')
#The time variable in the matrix refers to dicom volumes, and are found using the MRI repition time or TR of 1.5s
tr = 1.5
events_df['dicom_volume'] = (events_df['onset'] / tr).round().astype(int)
events_df.sort_values(by='onset')
| onset | duration | trial_type | response | dicom_volume | |
|---|---|---|---|---|---|
| 10 | 3.5 | 1.2 | banana | 2 | 2 |
| 11 | 32.3 | 1.4 | banana | 1 | 22 |
| 0 | 61.3 | 1.8 | orange | 4 | 41 |
| 5 | 89.8 | 1.6 | smelly-cheese | 1 | 60 |
| 16 | 121.9 | 1.2 | asafetida | 1 | 81 |
| 17 | 152.4 | 2.2 | asafetida | 1 | 102 |
| 1 | 180.2 | 1.2 | orange | 3 | 120 |
| 6 | 208.4 | 2.1 | smelly-cheese | 2 | 139 |
| 12 | 237.8 | 2.5 | banana | 2 | 159 |
| 7 | 264.7 | 2.9 | smelly-cheese | 2 | 176 |
| 18 | 294.0 | 2.5 | asafetida | 1 | 196 |
| 2 | 323.7 | 2.1 | orange | 4 | 216 |
| 13 | 353.3 | 2.7 | banana | 3 | 236 |
| 8 | 379.5 | 2.1 | smelly-cheese | 2 | 253 |
| 14 | 410.8 | 1.3 | banana | 3 | 274 |
| 19 | 437.5 | 1.1 | asafetida | 1 | 292 |
| 3 | 464.3 | 1.2 | orange | 4 | 310 |
| 4 | 491.8 | 1.3 | orange | 3 | 328 |
| 15 | 522.4 | 0.8 | banana | 2 | 348 |
| 9 | 551.7 | 1.6 | smelly-cheese | 2 | 368 |
# Define positive and negative response conditions
positive_conditions = ['orange', 'banana']
negative_conditions = ['isovaleric', 'bacon', 'asafetida', 'smelly-cheese']
# TR value
tr = 1.5
# Function to process events, filter by condition, and extract 4D data for 16 DICOM volumes after onset
def process_and_extract_data(file_path, nifti_file, condition_filter):
# Check if the NIfTI file exists
if not os.path.exists(nifti_file):
raise FileNotFoundError(f"NIfTI file not found: {nifti_file}")
# Extract subject and run from the file path
parts = file_path.split('/')
subject = parts[1]
filename = os.path.basename(file_path) # Get the file name from the path
run_modality = filename.split('_')[2].replace('_events.tsv', '') # Extract run_modality as 'run-01' or 'run-02'
# Ensure run_modality is in the format 'run-01' or 'run-02'
if run_modality not in ['run-01', 'run-02']:
raise ValueError(f"Unexpected run_modality: {run_modality}")
# Load the events.tsv file into a DataFrame
events_df = pd.read_csv(file_path, delimiter='\t')
# Check for 'trial_type' or 'condition' and rename to 'smell'
if 'trial_type' in events_df.columns:
events_df.rename(columns={'trial_type': 'smell'}, inplace=True)
elif 'condition' in events_df.columns:
events_df.rename(columns={'condition': 'smell'}, inplace=True)
else:
print(f"'trial_type' or 'condition' column not found in file: {file_path}. Adding 'smell' column with default value.")
events_df['smell'] = 'unknown' # Add 'smell' column with a default value if neither column is found
# Convert onset to DICOM volume
events_df['onset'] = events_df['onset'].astype(float)
events_df['dicom_volume'] = (events_df['onset'] / tr).round().astype(int)
# Filter based on condition
filtered_df = events_df[events_df['smell'].isin(condition_filter)].copy()
if filtered_df.empty:
return pd.DataFrame() # Return empty DataFrame if no valid events
# Load the NIfTI file for the run
img = nib.load(nifti_file)
img_data = img.get_fdata()
# Debugging: Check the shape of the NIfTI data
print(f"Loaded NIfTI file: {nifti_file}, Shape: {img_data.shape}")
# Get the shape of the data
_, _, _, num_volumes = img_data.shape
# Create a list to store the data for each event
extracted_data = []
for _, row in filtered_df.iterrows():
dicom_start = row['dicom_volume']
# Ensure we don't exceed the number of volumes in the run
dicom_end = dicom_start + 16
# Check if there are enough volumes to extract
if dicom_start >= 0 and dicom_end <= num_volumes:
# Extract the 4D data for the 16 volumes starting from the onset
data_slice = img_data[:, :, :, dicom_start:dicom_end]
# Only add to list if data_slice has at least 16 volumes
if data_slice.shape[-1] >= 16:
extracted_data.append({
'subject': subject,
'run_modality': run_modality,
'onset': row['onset'],
'dicom_start': dicom_start,
'dicom_end': dicom_end,
'smell': row['smell'], # Now using 'smell'
'data': data_slice
})
else:
print(f"Skipping event with fewer than 16 volumes: Start = {dicom_start}, End = {dicom_end}, Total Volumes = {num_volumes}")
else:
print(f"Skipping event with invalid DICOM volume range: Start = {dicom_start}, End = {dicom_end}, Total Volumes = {num_volumes}")
# Convert the list of dictionaries to a DataFrame
return pd.DataFrame(extracted_data)
# Process all subjects and runs for positive and negative responses
event_files = glob.glob('ds002185-download/sub-*/func/*_events.tsv')
positive_dfs = []
negative_dfs = []
for file_path in event_files:
# Assume the corresponding NIfTI file is named similarly to the events file
nifti_file = file_path.replace('_events.tsv', '_bold_preprocessed.nii.gz')
try:
positive_df = process_and_extract_data(file_path, nifti_file, condition_filter=positive_conditions)
negative_df = process_and_extract_data(file_path, nifti_file, condition_filter=negative_conditions)
# Append to the respective list if the DataFrame is not empty
if not positive_df.empty:
positive_dfs.append(positive_df)
if not negative_df.empty:
negative_dfs.append(negative_df)
except Exception as e:
print(f"Error processing file {file_path}: {e}")
# Combine all subjects into single DataFrames for positive and negative responses
if positive_dfs:
positive_combined_df = pd.concat(positive_dfs, ignore_index=True)
else:
positive_combined_df = pd.DataFrame()
if negative_dfs:
negative_combined_df = pd.concat(negative_dfs, ignore_index=True)
else:
negative_combined_df = pd.DataFrame()
# Print total number of events
print(f"Total number of positive events: {len(positive_combined_df)}")
print(f"Total number of negative events: {len(negative_combined_df)}")
Loaded NIfTI file: ds002185-download/sub-13/func/sub-13_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 405) Loaded NIfTI file: ds002185-download/sub-13/func/sub-13_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 405) Loaded NIfTI file: ds002185-download/sub-13/func/sub-13_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 414) Loaded NIfTI file: ds002185-download/sub-13/func/sub-13_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 414) Loaded NIfTI file: ds002185-download/sub-14/func/sub-14_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-14/func/sub-14_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-14/func/sub-14_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-14/func/sub-14_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-22/func/sub-22_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Loaded NIfTI file: ds002185-download/sub-22/func/sub-22_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Skipping event with invalid DICOM volume range: Start = 433, End = 449, Total Volumes = 440 Error processing file ds002185-download/sub-22/func/sub-22_task-odors_run-01_events.tsv: NIfTI file not found: ds002185-download/sub-22/func/sub-22_task-odors_run-01_bold_preprocessed.nii.gz Loaded NIfTI file: ds002185-download/sub-15/func/sub-15_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-15/func/sub-15_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-15/func/sub-15_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 434) Loaded NIfTI file: ds002185-download/sub-15/func/sub-15_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 434) Loaded NIfTI file: ds002185-download/sub-12/func/sub-12_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 460) Loaded NIfTI file: ds002185-download/sub-12/func/sub-12_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 460) Loaded NIfTI file: ds002185-download/sub-12/func/sub-12_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 460) Loaded NIfTI file: ds002185-download/sub-12/func/sub-12_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 460) Loaded NIfTI file: ds002185-download/sub-08/func/sub-08_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-08/func/sub-08_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-08/func/sub-08_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-08/func/sub-08_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-01/func/sub-01_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 415) Loaded NIfTI file: ds002185-download/sub-01/func/sub-01_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 415) Loaded NIfTI file: ds002185-download/sub-01/func/sub-01_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-01/func/sub-01_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-06/func/sub-06_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-06/func/sub-06_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-06/func/sub-06_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 429) Loaded NIfTI file: ds002185-download/sub-06/func/sub-06_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 429) Loaded NIfTI file: ds002185-download/sub-07/func/sub-07_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-07/func/sub-07_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-07/func/sub-07_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-07/func/sub-07_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-09/func/sub-09_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Loaded NIfTI file: ds002185-download/sub-09/func/sub-09_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Loaded NIfTI file: ds002185-download/sub-09/func/sub-09_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Loaded NIfTI file: ds002185-download/sub-09/func/sub-09_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 26, 440) Loaded NIfTI file: ds002185-download/sub-17/func/sub-17_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 425) Loaded NIfTI file: ds002185-download/sub-17/func/sub-17_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 425) Loaded NIfTI file: ds002185-download/sub-17/func/sub-17_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 425) Loaded NIfTI file: ds002185-download/sub-17/func/sub-17_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 425) Loaded NIfTI file: ds002185-download/sub-10/func/sub-10_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 432) Loaded NIfTI file: ds002185-download/sub-10/func/sub-10_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 432) Loaded NIfTI file: ds002185-download/sub-10/func/sub-10_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-10/func/sub-10_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-19/func/sub-19_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 417) Loaded NIfTI file: ds002185-download/sub-19/func/sub-19_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 417) Loaded NIfTI file: ds002185-download/sub-19/func/sub-19_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-19/func/sub-19_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-20/func/sub-20_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-20/func/sub-20_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-20/func/sub-20_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 421) Loaded NIfTI file: ds002185-download/sub-20/func/sub-20_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 421) Loaded NIfTI file: ds002185-download/sub-18/func/sub-18_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 421) Loaded NIfTI file: ds002185-download/sub-18/func/sub-18_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 421) Loaded NIfTI file: ds002185-download/sub-18/func/sub-18_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-18/func/sub-18_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-11/func/sub-11_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 428) Loaded NIfTI file: ds002185-download/sub-11/func/sub-11_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 428) Loaded NIfTI file: ds002185-download/sub-11/func/sub-11_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 413) Loaded NIfTI file: ds002185-download/sub-11/func/sub-11_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 413) Loaded NIfTI file: ds002185-download/sub-16/func/sub-16_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 463) Loaded NIfTI file: ds002185-download/sub-16/func/sub-16_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 463) Loaded NIfTI file: ds002185-download/sub-16/func/sub-16_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 465) Loaded NIfTI file: ds002185-download/sub-16/func/sub-16_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 465) Loaded NIfTI file: ds002185-download/sub-05/func/sub-05_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 414) Loaded NIfTI file: ds002185-download/sub-05/func/sub-05_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 414) Loaded NIfTI file: ds002185-download/sub-05/func/sub-05_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 420) Loaded NIfTI file: ds002185-download/sub-05/func/sub-05_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 420) Loaded NIfTI file: ds002185-download/sub-02/func/sub-02_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-02/func/sub-02_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-02/func/sub-02_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-02/func/sub-02_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-03/func/sub-03_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 398) Loaded NIfTI file: ds002185-download/sub-03/func/sub-03_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 398) Skipping event with invalid DICOM volume range: Start = 385, End = 401, Total Volumes = 398 Loaded NIfTI file: ds002185-download/sub-03/func/sub-03_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 404) Loaded NIfTI file: ds002185-download/sub-03/func/sub-03_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 404) Loaded NIfTI file: ds002185-download/sub-04/func/sub-04_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-04/func/sub-04_task-odors_run-01_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 440) Loaded NIfTI file: ds002185-download/sub-04/func/sub-04_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 422) Loaded NIfTI file: ds002185-download/sub-04/func/sub-04_task-odors_run-02_bold_preprocessed.nii.gz, Shape: (110, 104, 27, 422) Total number of positive events: 424 Total number of negative events: 410
# Display the DataFrames
print("Positive Events DataFrame:")
print(positive_combined_df.head())
print("\nNegative Events DataFrame:")
print(negative_combined_df.head())
# Count the number of events per subject in positive_combined_df
positive_event_counts = positive_combined_df['subject'].value_counts()
# Count the number of events per subject in negative_combined_df
negative_event_counts = negative_combined_df['subject'].value_counts()
# Print the znumber of events per subject for positive events
print("Number of positive events per subject:")
for subject, count in positive_event_counts.items():
print(f"Subject {subject}: {count} positive events")
# Print the number of events per subject for negative events
print("\nNumber of negative events per subject:")
for subject, count in negative_event_counts.items():
print(f"Subject {subject}: {count} negative events")
Positive Events DataFrame:
subject run_modality onset dicom_start dicom_end smell
0 sub-13 run-02 4.2 3 19 orange \
1 sub-13 run-02 153.0 102 118 orange
2 sub-13 run-02 269.7 180 196 orange
3 sub-13 run-02 383.7 256 272 orange
4 sub-13 run-02 507.5 338 354 orange
data
0 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
1 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
2 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
3 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
4 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
Negative Events DataFrame:
subject run_modality onset dicom_start dicom_end smell
0 sub-13 run-02 94.7 63 79 smelly-cheese \
1 sub-13 run-02 181.7 121 137 smelly-cheese
2 sub-13 run-02 241.7 161 177 smelly-cheese
3 sub-13 run-02 446.3 298 314 smelly-cheese
4 sub-13 run-02 568.4 379 395 smelly-cheese
data
0 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
1 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
2 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
3 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
4 [[[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ...
Number of positive events per subject:
Subject sub-16: 23 positive events
Subject sub-12: 22 positive events
Subject sub-17: 21 positive events
Subject sub-10: 21 positive events
Subject sub-03: 21 positive events
Subject sub-05: 21 positive events
Subject sub-11: 21 positive events
Subject sub-18: 21 positive events
Subject sub-20: 21 positive events
Subject sub-19: 21 positive events
Subject sub-04: 21 positive events
Subject sub-14: 21 positive events
Subject sub-07: 21 positive events
Subject sub-06: 21 positive events
Subject sub-01: 21 positive events
Subject sub-08: 21 positive events
Subject sub-02: 20 positive events
Subject sub-13: 20 positive events
Subject sub-15: 19 positive events
Subject sub-09: 18 positive events
Subject sub-22: 8 positive events
Number of negative events per subject:
Subject sub-16: 23 negative events
Subject sub-12: 22 negative events
Subject sub-13: 20 negative events
Subject sub-14: 20 negative events
Subject sub-02: 20 negative events
Subject sub-11: 20 negative events
Subject sub-18: 20 negative events
Subject sub-20: 20 negative events
Subject sub-19: 20 negative events
Subject sub-10: 20 negative events
Subject sub-17: 20 negative events
Subject sub-09: 20 negative events
Subject sub-07: 20 negative events
Subject sub-06: 20 negative events
Subject sub-01: 20 negative events
Subject sub-08: 20 negative events
Subject sub-15: 20 negative events
Subject sub-04: 20 negative events
Subject sub-05: 19 negative events
Subject sub-03: 19 negative events
Subject sub-22: 7 negative events
# Print the columns of the DataFrames to confirm
print("Positive DataFrame columns:", positive_combined_df.columns)
print("Negative DataFrame columns:", negative_combined_df.columns)
# Filter the DataFrames based on subject and run_modality
positive_filtered_df = positive_combined_df[(positive_combined_df['subject'] == 'sub-01') & (positive_combined_df['run_modality'] == 'run-01')]
negative_filtered_df = negative_combined_df[(negative_combined_df['subject'] == 'sub-02') & (negative_combined_df['run_modality'] == 'run-01')]
# Check if the filtered DataFrames are not empty and access the first row
if not positive_filtered_df.empty and not negative_filtered_df.empty:
# Get the data for the first positive and negative events
positive_data = positive_filtered_df['data'].iloc[0]
negative_data = negative_filtered_df['data'].iloc[0]
# Define the subject and run
subject = 'sub-01'
run = 'run-01'
# Define the NIfTI file path based on subject and run
nifti_file = f'ds002185-download/{subject}/func/{subject}_task-odors_{run}_bold.nii.gz'
# Check if the file exists
if os.path.exists(nifti_file):
# Load the NIfTI file
img = nib.load(nifti_file)
# Extract the affine transformation matrix
affine = img.affine
# Convert the extracted data to NIfTI images if needed
positive_img = nib.Nifti1Image(positive_data, affine)
negative_img = nib.Nifti1Image(negative_data, affine)
# Compute the mean image using nilearn
mean_positive_img = mean_img(positive_img)
mean_negative_img = mean_img(negative_img)
# View the mean image for an example negative
print("Viewing mean image for the selected negative event")
plotting.view_img(mean_negative_img, threshold=None, title='Mean Image for Selected Negative Event - Subject 01')
Positive DataFrame columns: Index(['subject', 'run_modality', 'onset', 'dicom_start', 'dicom_end', 'smell',
'data'],
dtype='object')
Negative DataFrame columns: Index(['subject', 'run_modality', 'onset', 'dicom_start', 'dicom_end', 'smell',
'data'],
dtype='object')
Viewing mean image for the selected negative event
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order)
# View the mean image for an example positive
print("Viewing mean image for the selected positive event")
plotting.view_img(mean_positive_img, threshold=None, title='Mean Image for Selected Positive Event - Subject 01')
Viewing mean image for the selected positive event
/usr/local/lib/python3.11/site-packages/numpy/core/fromnumeric.py:758: UserWarning: Warning: 'partition' will ignore the 'mask' of the MaskedArray. a.partition(kth, axis=axis, kind=kind, order=order)